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I,  INTRODUCTION 


The  emphasis  of  this  thesis  is  two-fold.  The  first  point 
is  to  gain  insight  into  the  characteristics  of  a spectral 
primitive  equation  model.  The  second  is  to  attempt  to  take 
advantage  of  the  spectral  expansions  in  an  initialization 
■ ec  hn i q ue . 

A spectral  model  uses  a Galerkin  formulation  with  the 
spherical  harmonics  as  its  choice  of  basis  functions.  The 
basis  functions  are  eigenfunctions  of  the  Laplacian  operator 
in  spherical  geometry  and,  thus,  are  used  to  advantage  for 
global  integrations.  An  immediate  consequence  of  the  Galerkin 
formulation  is  that  spatial  derivatives  of  a single  wave  are 
computed  without  the  usual  truncation  error  present  in  finite 
difference  methods.  This  means  that  the  spectral  model 
contains  no  linear  phase  speed  errors. 

The  non-linear  terms  are  computed  via  a transform  grid 
in  physical  space  (Orszag,  1971).  This  procedure  is  respon- 
sible for  making  the  efficiency  of  the  spectral  model  com- 
parable with  finite  difference  methods.  In  addition,  the 
non-linear  terms  may  be  computed  without  aliasing.  This  pre- 
vents non-linear  instability  and  conserves  invariant  integral 
properties.  Physical  processes  are  also  computed  on  the 
transform  grid. 

The  second  objective  of  this  thesis  is  to  test  a divergent 
initialization  for  the  spectral  model.  Most  static  and 
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variational  initialization  techniques  of  global  primitive 
equation  models  use  a constraint  to  balance  the  rotational 
part  of  the  wind  with  the  mass  field  while  setting  the 
divergent  part  to  zero.  Phillips  (I960)  has  shown  that  a 
completely  non-d i vergent  initialization  cannot  eliminate  the 
gravity  waves  which  would  be  generated  and  his  work  suggests 
that  the  quas i -geostroph i c divergence  could  tend  to  suppress 
these  gravity  wave  modes  for  mid-latitude  synoptic  scales. 

In  addition,  if  the  quas i-geostroph  i c divergence  had  internal 
modes  close  to  those  which  would  be  developed  by  t e model, 
such  things  as  the  large  scale  precipitation  would  be  more 
realistic  during  the  first  few  hours  of  the  model  forecast. 

It  is  an  aim  of  this  thesis  to  develop  a diagnostic 
velocity  divergence  to  insert  into  the  fields  of  an  otherwise 
completely  balanced  system  and  to  test  the  subsequent  effect 
on  the  spurious  gravity  wave  noise  and  the  developing  diver- 
gence. 

A spectral  formulation  of  the  primitive  equations  can 
lend  itself  to  techniques  which  initialize  according  to  scale. 
For  example,  the  quas i -geostroph i c divergence  might  be  com- 
puted only  for  the  synoptic  scales  of  the  model.  Furthermore, 
the  model  conveniently  supplies  diagnostic  information  about 
the  various  scales. 

The  final  goal  of  this  thesis  is  to  develop  a semi-impli- 
cit time  differencing  scheme  and  test  its  effects  during  the 
first  few  model  hours  after  initialization. 
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I I . MODEL  DESCR I PT ION 


The  model  described  here  has  been  developed  by  Rosmond 
(1977)  and  the  details  are  presented  here  for  reference 
purposes.  Similar  formulations  have  been  done  by  Hoskins 
and  Simmons  (1975)  and  by  Bourke  (1974).  The  equations  for 
an  inviscid,  adiabatic  and  hydrostatic  atmosphere  may  be 
written  as 

1.1)  = -7*(C+f )V  - k • Vx(RTVq  + a |X) 

1.2)  = k*Vx(^+f)V  - V*(RTVq  + - V ((f)  + j ) 

1.3)  . ve  - 

1.4)  |a= -D -V . Vq -|i 

1.5)  o||=-RT 

where 

C - vorticity  (^  = k • V x ^) 

D - veiocity  divergence  (D  = V • V) 

T - temperature 

6 - potent iai  temperature 

IT  - surface  pressure 

V - horizontal  velocity  vector 

(f)  - geopotential  height 

R - gas  constant 

Cp-  specific  heat 
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f - Coriolis  parameter 
o - vertical  coordinate  (o  = p/ir) 
a - vertical  velocity  (a  = do/dt) 
q - £n7T 

The  continuity  equation  (1.4)  may  be  refined  by  integrat- 
ing in  the  vertical  and  imposing  boundary  conditions  on  o; 

• • 

0(0)  = 0(1)  = 0.  Introducing  the  notation  G = V*Vq,  we  may 
rewr ite  Eq.  (1.4)  as 

1.6)  -^  = -(D  + G)  . 

where  the  over  bar  denotes  a vertical  average. 

The  vertical  velocity,  o,  may  be  obtained  diagnostically 
by  substituting  Eq.  (1.6)  into  Eq.  (1.4)  and  integrating  in 
the  vertical  to  obtain 


-/ 


.7)  o = (D  + G)o  - / (G  + D)  do  . 

o 

By  defining  a horizontally  mean  temperature, 


T = T (o)  + T'(o,  X,  6,  t) 


and  using  the  following  operator 


/V  I 3X  ^ 3Y 

I-P 
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We  can  now  expand  the  vector  quantities  in  our 
equations  in  spherical  coordinates  (X,n)  to  arrive 


.8) 

3C  _ 
at 

-a(A,B) 

.9) 

3D 

at 

a(B,-A)  -V^Cr^E 

« 

+ ij)  + RT  q ) 

10) 

ae 

at 

-a(Ue,V0)  + 0D  - 

. 30 

° 30 

1 1 ) 

II 

-(D  + G) 

12) 

3d)  _ 
° 30 

-RT 

where 


A 

B 

G 

E 

U 

V 

X 

n 

u 


• 8V 

= (c+f)u  + + 


RT'  „ 
-y  cos  0 

r 


= (c+f)v  - 6 


U dZnn 

2 3X 

-U 


7~ 

2( l-M  ) 

= u cos  n/r 


^ _ RT' 

3o  ^ 


+ V 


dZnir 

3y 


3S,nTr 

~TT 


= V cos  n/r 


= longitude 
= latitude 


= si n n 


basic 

at 
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Equations  (1.8)  - (1.12)  are  the  basic  equations  used  in 
the  model.  The  equations  are  represented  spectrally  in  the 
horizontal  and  finite  differenced  In  the  vertical.  The  de- 
pendent variables  are  written  in  terms  of  a triangularly 
truncated  series  of  spherical  harmonics; 

X = z xj  yJ 

where  the  summation  is  a double  sum  over  m,n  for  |m|  _<  M and 
|m|  £ 2.  £ M.  The  separation  is  such  that  the  coefficients, 
X^,  are  functions  of  time  and  vertical  coordinate  and  the 
spherical  harmonics,  Y^,  are  horizontal  functions  of  space. 
The  normalization  and  orthogonality  properties  of  the 
allow  the  coefficients  to  be  obtained  as  follows: 

2it  +1 

= ‘X'  ''J>  = i / / X x;  dMdX  . 

o - 1 

The  model  used  for  this  study  has  five  vertical  layers 
staggered  as  in  Fig.  I.  The  even  levels  carry  o and  the  odd 
levels  carry  all  the  dependent  variables  plus  the  diagnostic 
variables,  U,V. 

The  finite  difference  scheme  for  the  vertical  advection 
terms  is 

Y I • • 

0 = TTir — r 1^0.]. 

9a  2(Aa,  ) k+ 1 k k+l  k k- 1 k 

k 

This  formulation  conserves  kinetic  energy. 


15 


I.D.T.U.V 

Ld.t.u.v 

» 

^2. 

<iaL!JJ£. 

64 

03 

^E£,LyjL 

^5 

VERTICAL  SCHEME 
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Following  Arakawa  and  Lamb  (I976),  the  hydrostatic  equa- 
tion (I.I2)  can  be  finite  differenced  as 

^k+\  ■ '<’k  = ■ ^ Vi  ®k  "^k’  j 

where 

“k  = ' \ - ' • 

A boundary  condition  is  formed  by  integrating  the  hydrostatic 
equation  from  a = 0 to  I which  gives 


o o 


where  is  the  -terrain  geopotential.  This  boundary  con- 

dition is  then  finite  differenced  as 

1.14)  I *^(40^)  - = R £ 

Equations  (I.I3)  and  (1.I4)  are  then  combined  to  give  the 
matrix  equation 

1. 1 5)  (i)  = [c:t  + 

which  is  the  finite  difference  form  of  the  hydrostatic  equa- 
tion where  the  variables  are  now  column  vectors. 
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The  pressure  tendency  equation  (l.ll)  can  be  written  as 


|a  . - E(G^  ♦ 0,)(A.,,> 

k 

or  in  matrix  form  as 

1.16)  = -[N]  (G+D) 

Similarly  equation  (1.7)  can  be  written  in  matrix  form  as 

1.17)  a = [Z]  (G+D) 

The  thermodynamic  equation  (I.IO)  is  differenced  consis- 
tent with  the  method  employed  in  the  hydrostatic  equation. 

The  potential  temperature  is  written  in  terms  of  temperature 
and  the  last  term  becomes 


- ■'k>* 


“k.i  Vi 


(G.  - G - 0). 
k 


For  the  purposes  of  the  sem i - i mp I i c i t formulation  de- 
scribed later  in  the  paper,  we  separate  the  temperature  into 
its  horizontal  average,  T , and  perturbation,  T*.  The 
matrix,  [yJ,  is  defined  such  that 
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The  vertically  averaged  divergence  can  be  written  6 
[N]D,  so  that 

-CY]a  - ^ D = -CYDi  - ^ [nDD  . 
p p 

Substituting  for  o from  equation  (1.17)  we  obtain 
* » 

-CYDa  - ^ [N]D  = -CY]CZ]G  -([Y][z:  + gL  [n])D 
P P 

= -[s]G  + [q:d  . 

Finally  the  entire  term  can  be  written  as 


.18)  6 “ . -CS]G+C«D-  - T') 

Vl 

o,  R/C  RT,  RT' 

* c^‘ V®’-  ^ 5 

“k-l  p " 


The  nonlinear  terms  are  computed  using  the  transform 
method  suggested  by  Orszag  (1971).  The  longitudinal  direc- 
tion is  done  using  Fast  Fourier  Transform  and  the  latitudinal 
direction  is  done  with  Gaussian  Quadrature.  The  number  of 
latitudes,  N,  and  longitudes,  M,  satisfy: 

N _>  (3J  + I )/2  and  M ^ 3J  + I . 

The  number  of  points  are  picked  to  ensure  non-aliased  results 
for  quadratic  terms.  Certain  terms  containing  a have  triple 
products  and  are  not  computed  alias  free.  Experiments  by 
the  author  and  a previous  study  (Hoskins  and  Simmons,  1975) 
indicate  that  this  source  of  error  is  negligible  for  the 
length  of  the  integrations  performed  in  this  study. 
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III.  SEMI-IMPLICIT  SCHEME 

This  section  gives  details  of  the  sem i - i mp I i c i t time 
differencing  scheme  used  in  some  of  the  experiments  of  this 
study.  Using  the  matrices  defined  in  the  model  description 
section,  equations  (1.9)-(1.I2)  can  be  written 


3D 

TF 


2.1)  ^ + RT  q)  = Fr 


2.2)  - [Q]D  = 


2.3)  + [N]D  = -G 


2.4)  (}.'  = CCjT 


where 


^ = 0 — (p 


s f c 


Fp  and  F.^  represent  the  remaining  terms  in  each  equation 
which  are  not  explicitly  separated  out.  Following  Robert 
et  al.  (1972),  all  the  terms  on  the  left  hand  side  of  equa- 
tions (2.l)-(2.4)  are  time  averaged.  We  define  the  following 
averaging  operator  and  represent  the  local  tendencies  by  a 
Leapfrog  finite  difference  scheme: 

( ) ; 


TF  ‘ > - > ' 


( 

2TSt) 
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Equations  (2.1)-(2.4)  can  then  be  written  as 


2.5) 

6^D  = 

Fj^  - 7 (^  + RT 

2.6) 

+ CQ] 

2.7) 

6^q  = 

- G - CN] 

2.8) 

= 

Cc] 

El  i m i n a t i n g in  equation  (2,5)  by  substitution  of  Eq. 

(2.8),  we  arrive  at  a time  averaged  set  of  equations: 

+ (At)(F|^  - V^([Cjr^  + RT*  q^)) 

+ (At)(F^  + [Q]  d'^) 

-(At)(G  + [N]  d"^)  . 


2.9) 

d'^ 

= 

2.10) 

jT-At 

2.11) 

t-At 

= q 

Substituting  equations  (2.11)  and  (2.10)  into  (2.9)  and 
solving  for  we  obtain  a Helmholtz  equation  for  the  time 
averaged  divergence 

2.12)  ([l>(At)^[B]7^)D^  = + (At)[Fp-7^([C]T‘^‘^''’ 

+ (At)[C]F^+RT*q‘'’"'^''’-(At)  RT*  G] 

where  Cl]  is  the  identify  matrix  and  Cb]  = [C][0]  - RT  [N]. 

For  a spectral  formulation,  the  equation  is  trivial  to 

solve  since  the  spherical  harmonics  are  eigenfunctions  of 
2 -t 

the  Eaplacian  7 . Knowing  D , we  then  compute  the  prognostic 
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variables  from  equations  (2.6)  - (2.7). 

is  computed  from  the  definition  of  the  time  averaging 


operator : 


. 2D+  - 


and  the  vorticity,  never  entered  the  sem  i - i mp  I i c i t 

formalism  and  is  computed  from  a standard  Leapfrog  scheme. 
Calculations  with  the  sem i - i mp I i c i t model  required  5$  more 
computer  time  per  time  step.  Stable  I nteg>-at  1 ons  were 
possible  with  a time  step  of  60  mins. 


1 V . QUAS I-GEOSTROPHIC  DIVERGENCE 


Following  the  analysis  of  Simmons  and  Hoskins  (1976),  a 
diagnostic  divergence  Is  obtained  for  the  case  of  quasi- 
geostrophic  flow  with  the  full  variation  of  the  Coriolis 
parameter  over  the  sphere.  The  relevant  equations  become 

-It'  - ^ 


3.2)  JC;  - V^(4>'  + RT*q)  = 0 


3.3,  - = 1 -|l‘, 


where 


V 

0 

/D 


rotational  part  of  V 
sigma  dot  at  odd  levels 
V • (fV(V“^D)) 


Making  use  of  the  hydrostatic  relationship,  equation 
(1.15),  we  may  write  equation  (3.3)  as 

If  = -CC](V^  • VT)  + Y*  i 
* # 

where  y*  = I - ). 

2 3 

Taking  V of  this  equation,  if  of  equation  (3.1)  and 

of  eq.  (3.2)  and  substituting  we  arrive  at  the  following 

A 

equat  i on  for  o: 

3.4)  Y*o  = -v'^C  (V^*V(i;+f))+  ^D]+CC](V^*VT)+RtV  . 


A second  equation  for  o is  obtained  from  the  continuity 


equation  in  the  form  of  eq,  (1.17),  which  when  explicitly 
wr i tten  out.  Is 


3.5)  0 


0(5  + 5) 


a 

(D  + G)  do 


o 

where  o Is  specified  on  even  levels. 

/V 

The  diagnostic  o from  equation  (3.4)  can  be  used  to  com- 
pute the  divergence  using  eq.  (3.5)  provided  we  specify  a 
boundary  condition.  For  the  purposes  of  Initialization,  this 
is  most  conveniently  chosen  to  be 


1^  = - (5  + G)  = 0 


since  this  will  tend  to  eliminate  external  gravity  waves. 
Applying  this  constraint  to  equations  (3.4)  and  (3.5)  we  ob- 
tain 

3.6)  = -V“^(  (V^^*V(C+f))+  <^D)  + CCDV^'VT 

/o 

(D  + G)  do  . 
o 

A 

The  method  of  solution  will  be  to  compute  o using  equa- 
tion (3.6)  from  an  initial  guess  of  the  divergence.  The  even 
level  o’s  are  then  found  by  interpolation  and  a new  guess  of 
the  divergence  is  obtained  by  inverting  equation  (3.7).  The 
updated  divergence  can  be  substituted  back  into  equation  (3.6) 
and  the  procedure  continues  iteratively  until  sufficient 
accuracy  is  obtained. 
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Equations  (3.6)  and  (3.7)  will  now  be  put  Into  spectra 
form.  The  equations  contain  three  similar  advective  terms, 
It  will  be  sufticlent  to  show  the  spectral  form  of  one  of 


them : 


5.8)  V^-5(C.f)  . +t)  . -ij  . 


First,  the  terms  are  Individually  transformed  to  physica 
space  by  the  summations: 


dip  . m ° 0 1 mA 

£,m 

£,m 

dip  . , m nf'f'  n I 

S=f, i=t  ‘ 


3(^+f ) 

3y 


3P 

^^m^^ImA  ^ . 

£7m 


The  multiplications  In  eq . (3.8)  are  then  computed  at 

grid  points  and  the  result  is  transformed  back  to  spectral 

space  by  computing  the  inner  product  defined  as: 

+ 1 -2tt 


3.9)  = i 

L 4tt 


// 


J(ip,  C+f)vJ'  dAdy  = <J(i|;,  C+f),  Y^>  . 


The  integral  is  computed  using  Gaussian  quadrature  in 
the  y direction  and  Fast  Fourier  Transforms  in  the  A direc- 
tion. The  first  term  of  equation  (3.6)  can  then  be 
eva I uated : 

»C(V^‘V(i;+f))  = V *(fV(V^J)) 
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and  expanded  in  spectral  form  using  the  previous  computation 
of  J^: 

/'  +1  ■ 9P^ 

< (Vi|^*V(^+f))  yJ^>  = / UILJ  P^)PLdy-fi y* (|V)(I  ' 

where  the  X integration  has  already  been  done. 

Making  use  of  the  recurrence  relationship 


upM  = eM  M ^ M M 


" 41  - 


we  can  write  the  first  integrand  as 


j.  . ,M  pM  _ M M M M M M 

^ ^ ^ '^£-1  ^ '^£+1  ^£t|  ^£ 


and  using  the  orthogonality  of  the  Legendre  polynomials  we 
obtain  for  the  first  integral 
/' 

, / /r-  ^ ,M  M ^ ,M  M 

° J ^ ’ 

-I 

The  second  integral  can  be  done  similarly  using  the 
relationship 


3P^ 

^£  ^£-1  ■ ®'^£t|  ^£+1 
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to  obtain 


+ 1 

5.  I I )y*  ( l-y 


■M  -pM  M M 


Xn+TT  3y  L 


L+l  L+l  "L-l  L 


Combining  eqs  (3.10)  and  (3.11)  we  obtain 


1 iT\  nr  L+l  M L M ,M  -i 

3.12)  <jr(v^-V(ctf).Y|_>  . (!C-^e^  J|__l  J^^|]  . 


The  second  order  operator  tr^  can  now  be  formulated. 
Using  the  notation,  <jCd,  Y^>  = , we  can  write  by  analogy 

with  equation  (3.12) 


pM  _ pL+l  , L M -1 

^ Dl-1  ^i+r^L+i  '^L+i^  • 


Thus, 


<D  = E 
£ , m 


and 


<<^o.  y">  = <«»D>.  r">  - acii  F«_,  F ^ . 


Substituting  in  for  F^_j  and  F^^|,  we  obtain  the  result 


T 173  „2r,(L+l)L  M M 


+ ( 


+ ( 


L^-l  .,M.2  , L^+  2L  , M .2..M 

L (L+l)  ^ 

+ L M M M 
(L+KL+2)  ^L+l  ^L+2  °L+2-^  ’ 
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V 


DESIGN  OF  THE  EXPERIMENTS 


The  initial  conditions  used  in  the  experiments  of  this 

thesis  were  developed  from  a ba roc  I i n i ca 1 1 y unstable  mean 

flow  and  a small  perturbation  in  zonal  wave  no.  6 (m=6). 

The  mass  field  of  the  mean  flow  was  constrained  to  be  in 

geostrophic  thermal  balance  with  the  mean  wind.  The  vertical 

profile  of  the  mean  wind  was  linear  varying  from  4.5  m/sec 

at  o = I to  49.5  m/sec  at  a = 0.  The  latitudinal  structure 
2 

varied  as  sin  (2f|)  placing  a jet  max  at  n = 45®.  The  verti- 
cal temperature  profile  approximated  the  U.  S.  Standard 
Atmosphere . 

The  perturbation  in  wave  no.  6 corresponded  to  a maximum 
'.-component  of  .1  m/sec.  The  model  was  integrated  for  192 
hrs  allowing  the  baroclinic  wave  to  grow  to  finite  amplitude. 
Fig.  2 shows  a time  series  of  the  vorticity  during  this 
period.  The  growth  is  strikingly  linear  until  192  hrs  where 
there  is  evidence  of  the  mean  flow  (m  = 0)  changing  due  to 
feedback  from  the  eddies.  The  I inear  growth  corresponds  to 
a doubling  time  of  roughly  26  hrs  for  m = 6.  The  fields  at 
192  hrs  provided  the  balanced  fields  for  the  initialization 
exp er  i men t s . 

Three  types  of  experiments  and  one  control  case  were 
integrated  for  72  hrs  of  model  time. 
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Figure  2. 


Vo rt i c i ty  time 
Vertical  axis 
Solid  c u rve  i s 
latitude  60°N. 
dashes  are  m=2 


series  during  the  first  192  hrs. 
s logarithmic.  Units  are  secs“l 
fhe  magnitude  of  wave  no.  0 at 
Dotted  curve  is  m=6,  short 
and  long  dashes  are  m=l8. 
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The  control  case  consisted  of  initializing  with  the 
balanced  fields  and  integrating.  The  experiments  consisted 
of  the  following: 

i Experiment  A:  Initialize  with  the  balanced  fields  but 

; nodi vergence 

Experiment  B:  Initialize  with  the  balanced  fields  plus 

; the  quas i -geostroph i c divergence 

Experiment  C:  Identical  to  Experiment  A buf  with  a semi- 

imp I i c i t scheme 

I All  experiments  and  the  control  were  first  integrated 

1 with  a spectral  truncafion  of  wave  no.  21.  These  will  be 

referred  to  as  low  resolution  experiments.  High  resolution 
experiments  were  integrated  for  the  control  and  Experiments 
A and  B.  These  cases  were  truncated  at  wave  no.  42.  Time 
' steps  of  12  mins  for  the  low  resolution  and  6 mins  for  the 

high  were  used.  Experiment  C was  integrated  with  time  steps 
of  12  and  60  mins.  All  experiments  had  a weak  Robert  time 
filter  of  0.1  and  no  dissipation  term  was  used. 


31 


VI,  RESULTS 


A.  CONTROL  CASE 

Figures  3-6  show  the  initial  fields  at  level  5 for  the 
low  resolution  experiments.  Figures  7-10  are  the  correspond- 
ing high  resolution  fields.  They  are  characterized  by  a 
developing  frontal  zone.  The  advantage  of  the  high  resolu- 
tion is  apparent  in  the  temperature  gradients  across  the 
cold  front  and  warm  front  and  in  the  surface  pressure.  The 
low  in  the  high  resolution  surface  pressure  field  is  5 mb 
deeper.  The  negative  and  positive  areas  of  the  divergence 
and  vorticity  fields  are  approximately  symmetrical  in  the 
low  resolution  fields  and  quite  asymmetric  in  the  high  reso- 
lution fields.  The  region  of  maximum  positive  vorticity 
corresponds  very  well  with  the  maximum  temperature  gradient. 

The  mean  fields  consist  of  an  indirect  cell  induced  by 
the  growing  baroclinic  eddies.  The  indirect  cell  develops 
in  order  to  decrease  the  vertical  wind  shear  which  must  remain 
in  approximate  thermal  balance  with  the  decreasing  pole  to 
equator  temperature  gradient.  The  effect  is  more  pronounced 
in  the  low  resolution  fields.  The  high  resolution  mean 
fields  contain  a weak  indirect  circulation  and  a strong 
d i rect  cell. 

The  dynamic  situation  during  the  period  of  the  control 
case  is  one  of  a rapid  cascade  of  energy  to  smaller  scales. 
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Figure  8.  Initial  high  resolution  temperature 
interval  is  5 ° K. 
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DIVERGENCE  IN  1/SECS 


The  horizontal  scale  of  the  frontal  zone  decreases.  Figures 
11-13  show  this  effect  in  the  high  resolution  temptrature 
fields  for  times  out  to  36  hr.  Beyond  this  point,  the 
energy  Puilds  up  at  the  spectral  limit  of  the  model  since 
there  is  no  dissipation  and  the  72  hr  fields  are  non-meteoro- 
logical  (Fig.  14  and  15).  Because  the  energy  blockage  occurs 
at  smaller  scales,  the  high  resolution  case  (Fig.  14)  is, 
paradoxically,  noisier  and  in  greater  error. 

B.  GRAVITY  WAVES 

In  anticipation  of  later  sections,  we  discuss  the 
characteristics  of  inertial  gravity  waves  for  the  model  used 
in  this  study.  The  divergence  is  the  most  sensitive  indica- 
Tor  of  gravity  wave  motion.  As  discussed  in  Hoskins  and 
Simmons  (1975),  equation  (2.12)  is  the  sem i - i mp I i c i t analogue 
of  the  gravity  wave  equation  for  a multi-layer  model.  The 
eigenvalues  of  the  matrix  CB]  (eq.  2.12)  give  the  permissible 
e i gen ve I oc i t i es  of  the  gravity  wave  modes.  The  corresponding 
eigenvectors  give  the  vertical  structure.  Figure  16  depicts 
the  eigenmodes  for  the  experiments  in  this  study.  Table  I 
gives  .he  periods  of  the  gravity  wave  modes  for  zonal  wave 
no . 6 in  h r s . 


TABLE  I 
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VERTICAL  DIVERGENCE  MODES 


The  hotel  divergence  may  be  thought  of  as  consisting  of 
a sum  of  the  individual  eigenmodes: 

4.1)  D = I 

k 

where  the  D are  the  e i gen vecto rs . The  weighting  coefficients, 
a^^,  may  be  obtained  from  the  simultaneous  equations  formed 
by  the  dot  products  of  the  individual  eigenvectors  with  equa- 
tion (4.1).  A time  series  of  the  coefficients  will  provide 
useful  diagnostic  information  in  later  sections.  Figures  17 
and  18  are  a time  series  of  the  external  and  first  internal 
modes  for  the  low  resolution  control  case  at  a point  on  the 
60°  latitude  circle.  It  is  apparent  that  the  coefficients 
contain  both  the  meteorological  modes  and  any  excited  gravity 
modes.  Figure  17  contains  no  high  frequency  oscillation  and 
is  purely  meteorological.  Figure  18  contains  a meteorologi- 
cal mode  with  a period  of  roughly  3 days  plus  a smaller 
amplitude  wave  whose  period  of  approximately  1/2  day  corre- 
sponds well  with  the  calculated  value  for  the  first  internal 
mode  from  Table  I.  The  remaining  internal  gravity  modes 
have  periods  approaching  meteorological  values  and  are  diffi- 
cult to  distinguish. 

C.  EXPERIMENTS  A AND  B 
I . Low  Resolution 

Neglecting  the  divergence  in  the  initialized  fields 
(Experiment  A)  generated  comp u ta t i ona  I gravity  waves.  Compar- 
ing the  pressure  tendencies  of  the  control  case  (Fig.  I9c) 
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Figure  17.  Divergence  time  series  of  the  external 
mode  for  the  low  resolution  control, 
case.  Vertical  axis  units  are  sec”  . 
Time  is  In  hrs. 


with  those  of  Experiment  A in  Fig.  19a,  it  is  apparent  that 
high  frequency  oscillations  of  roughly  3 hr  period  are 
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excited.  The  initial  magnitudes  of  the  oscillations  are 
approximately  1-2  mb/hr.  The  oscillations  dampen  out  with 
time  as  the  geostrophic  adjustment  takes  place  although  the 
dispersal  of  the  gravity  waves  during  the  adjustment  process 
is  hampered  by  the  longitudinal  periodicity  and  the  symmetry 
of  the  fields  about  the  equator. 

Figure  19b  shows  the  pressure  Tendencies  with  the 
quas i -geostroph i c divergence  inserted  into  the  initializa- 
tion fields  (Experiment  B).  The  highest  frequency  oscilla- 
tions have  been  reduced  over  Experiment  A by  a factor  of 
a p p rox i mate  I y 2 . 

Comparison  of  the  divergence  modes  in  Figs  20,  21,  22 
allows  us  to  be  more  quantitative.  The  figures  have  the 
control  case  subtracted  out  and  contain  the  gravity  wave 
noise  generated  by  the  approximations  of  each  experiment. 

The  noise  generated  in  the  internal  modes  by  the  non-d i vergent 
initialization  is  an  order  of  magnitude  larger  'han  the  ex- 
ternal mode  and  is  generally  2b%  of  the  total  divergence. 

The  figures  show  that  the  external  mode  of  Experiment  B is 
about  68?  of  Experiment  A.  The  first  internal  mode  shows 
more  noise  in  Experiment  B initially  and  then  decreasing  to 
50?  and  the  second  internal  mode  shows  a decrease  of  about 
50?  over  Experiment  A.  In  general  the  effect  of  the  quasi- 
geostrophic  divergence  is  to  reduce  the  gravity  wave  noise. 


il 
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Figure  21.  Gravity  wave  noise  in  the  first  internal  mode  - otherwise 


2 . High  Resolution 

Figure  23  shows  the  quas i -geostroph i c divergence  at 
level  5 which  was  used  for  Experiment  B.  A qualitative  com- 
parison with  the  actual  model  divergence  (Fig.  10)  shows  the 
phases  to  be  in  agreement  but  that  the  q ua s i -geost rop h i c 
approximation  is  an  overestimate  for  the  perturbation  fields. 
The  zonal  mean  fields  do  not  compare  well.  The  quasi- 
geostrophic  divergence  shows  a strong  indirect  circulation 
which  is  displaced  poleward  relative  to  the  model  divergence. 
There  is  no  evidence  of  a direct  cell. 

Figures  24-26  depict  the  noise  generated  in  the 
divergence  modes  for  Experiments  A and  B.  General ly  they 
show  that  the  gravity  wave  noise  of  Experiment  B is  66^-72? 
of  that  in  Experiment  A.  The  effect  of  the  quasi-geostrophic 
divergence  is  to  reduce  the  noise.  However,  this  reduction 
is  less  than  in  the  low  resolution  experiments. 

In  view  of  the  fact  that  the  high  resolution  experi- 
ments were  less  successful  than  the  low  resolution  and  that 
quasi-geostrophic  formulations  are  valid  only  for  small 
values  of  the  Rossby  number,  it  is  interesting  to  compare  the 
initial  model  divergence  of  the  control  case  with  the  quasi- 
geostrophic  divergence  in  a spectral  manner.  Referring  to 
Fig.  27,  it  is  seen  that  there  is  good  agreement  between  the 
model  and  quasi-geostrophic  divergences  for  zonal  waves 
corresponding  to  m = 6,  12,  18.  The  overestimation  by  the 

quas i -geostroph i c divergence  noted  earlier  is  primarily  in 
wave  no.  6.  The  mean  fields  (m=0)  and  waves  corresponding 


to  m 


24,30,36,42  do  not  compare  well. 
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igure  26.  Gravity  wave  noise  in  the  second  interval  - otherwise 
as  in  Fig.  24 . 


One  motive  for  initializing  with  the  q u a s i -geo st ro- 
phic  divergence  is  tiiat,  presumably,  the  "true"  divergence 
would  be  approximated  at  the  start  and  that  the  model  would 
thus  develop  this  divergence  field  earlier  than  a non- 
divergent  initialization.  Figures  28  and  29  address  this 
point.  They  show  that  at  6 hours  into  the  integrations,  zonal 
wave  numbers  6,  12,  18  do,  in  fact,  approach  the  control 

divergence  slightly  faster.  There  is  almost  no  difference 
for  the  higher  wave  numbers  and  the  mean  field  shows  poor 
adjustment  in  either  case.  This  might  be  expected  for  the 
mean  field  since  the  geostrophic  adjustment  process  is  very 
inefficient  at  this  scale  in  a bounded  region.  On  the  other 
hand,  the  higher  wave  numbers  adjust  very  rapidly  and  after 
6 hrs,  they  are  in  good  agreement  with  the  control  case. 

In  light  of  the  results  showing  good  agreement  in 
the  initial  quas i -geostroph i c divergence  with  the  control 
divergence  for  zonal  wave  number  6,  12,  18  and  bad  agreement 

for  other  scales,  it  is  natural  to  modify  Experiment  B such 
that  only  scales  corresponding  to  m = 6,  12,  18  are  used 

in  the  qua 5 i -geost rop h i c divergence.  Figure  30  shows  the 
time  series  of  the  external  divergence  mode  for  this  experi- 
ment. Figure  31  shows  the  spectral  adjustment  of  the  diver- 
gence at  6 hours.  The  corresponding  results  for  the  unmodi- 
fied Experiment  B are  Figs.  24  and  29.  The  results  are 
almost  identical  to  the  unmodified  quas i -geostroph i c diver- 
gence results. 

( 
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Figure  30. 


T i me  series 
in  the  mod  i f 
scale  is  sec 


of  the  gravity  wave  noise 
i ed  Experiment  B.  Vertical 
Time  is  in  hrs. 


D. 


EXPERIMENT  C 


As  discussed  in  Hoskins  and  Simmons  (1975),  a linear  sta- 
bility analysis  of  the  sem i - i mp I i c i t scheme  shows  that  the 
time  step  is  limited  only  by  the  Rossby  waves.  The  effect 
of  the  scheme  is  to  slow  down  the  faster  moving  gravity  wave 
modes  so  that  the  integrations  do  not  amplify  in  time. 

Slower  moving  atmospheric  motions  are  less  affected.  Table  I 
presents  sample  calculations  of  the  distortion  of  the  periods 
of  the  gravity  waves  present  in  these  experiments  for  two 
time  steps.  These  may  be  compared  with  Table  I.  The  general 
behavior  is  that  the  periods  increase  with  increasing  time 
steps. 


TABLE  I I 


At 

Externa  1 

First 

Second 

Third 

Fourth 

12  mins 

3.31 

11.14 

28.0 

6 1.7 

154.0 

60  mins 

5.7 

12.2 

28.5 

6 1.9 

154.1 

Two  points  will  be  investigated  in  this  experiment.  The 
first  is  that  the  sem  i - i mp I i c i t scheme  should  tend  to  smooth 
out  high  frequency  oscillations,  particularly  with  a large 
time  step.  The  second  is  that  the  sem i - i mp I i c i t method 
might  possibly  hamper  the  geostrophic  adjustment  process 
since  it  slows  down  the  gravity  waves  which  are  the  mechanism 
for  the  adjustment. 
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Figure  32  shows  the  pressure  tendencies  using  the  semi- 
implicit  method  with  +lme  steps  of  12  and  60  mins.  Compar- 
ing this  with  the  corresponding  results  for  the  explicit 
integrations  (Fig.  19a),  it  is  evident  that  the  ef-^ect  of 
the  sem i - i mp I i c i t method  is  to  smooth  out  the  3 hr  oscilla- 
tions. Figure  33  shows  the  relative  effect  of  a sma  I I vs. 
large  time  step  in  the  sem i - i mp I 1 c i t method  on  the  external 
gravity  waves.  The  larger  period  of  the  60  min  time  step  is 
evident.  The  larger  time  step  initially  amplifies  the  com- 
putational noise,  but  after  24  hrs  it  is  considerably  damped. 
The  smoothing  in  the  pressure  tendencies  is  postulated  to 
come  about  from  the  decreased  velocity  of  the  external 
gravity  mode  assuming  an  advective  time  scale. 

The  second  point  of  this  experiment  is  shown  in  Figs.  34 
an^  35.  Comparison  of  the  divergence  coefficients  for  time 
steps  of  12  and  60  mins  shows  no  h'fference  at  3 hrs.  The 
adjustment  of  the  divergence  to  the  control  case  is  apparent- 
ly the  same  for  either  time  step. 
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vjiavi.y  Wdve  noise  in  rne  externa;  mode  for  Experimen 
with  time  steps  of  12  mins  (33a)  and  60  mins  (33b). 
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Figure  35 . 


3 hr  divergence  coefficients  with  a 60 
min  time  step  - otherwise  as  in  Fig,  34. 
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VII.  SUMMARY  AND  CONCLUSIONS 

The  results  of  this  study  are  based  upon  the  divergence 
developed  by  a baroclinic  wave  on  a uniform  earth.  Initiali- 
zation without  the  divergent  wind  component  generated  spurious 
inertial  gravity  waves.  The  internal  state  of  the  initial 

fields  had  nearly  zero  mean  divergence.  Thus,  the  external  i 

mode  of  the  generated  gravity  waves  was  an  order  of  magni-  I 

tude  less  than  the  total  divergence  and  was  not  significant.  | 

However,  the  internal  gravity  wave  modes  generated  were  signi-  j 

ficant.  An  important  source  for  the  external  mode,  namely  1 

topography,  was  neglected. 

Initializing  with  the  q ua s i -geost rop h i c divergence  | 

generally  decreased  the  spurious  gravity  wave  noise  and  aided  : 

the  ensuing  integrations  to  develop  a divergence  closer  to 
the  control  case,  but  both  effects  were  small.  Furthermore, 
the  advantage  of  the  q ua s i -g eost rop h i c divergence  decreased 
when  smaller  scales  were  included  in  the  integrations.  In- 
serting the  q ua s i -geo s t rop h i c divergence  into  only  synoptic 
scales  gave  no  advantage.  However,  the  effect  of  this  experi- 
ment may  have  been  reduced  by  the  fact  that  there  were  no 
planetary  scale  waves. 

Using  a sem i - i mp I i c i t scheme  with  the  same  time  step  as 
an  explicit  integration  had  no  effect  on  the  divergence  modes  | 

but  did  smooth  out  oscillations  in  the  pressure  tendencies.  ' 

i 

i 

i 

i 
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A larger  time  step  Initially  increased  the  amplitude  and 
period  of  the  fast  external  gravity  waves  and  had  less  effect 
on  slower  moving  motions  while  still  stabilizing  pressure 
tendencies. 

In  general,  the  effect  of  a divergent  initialization  for 
a global  primitive  equation  model  is  small.  This  same  con- 
clusion has  been  reached  by  two  similar  studies  - Houghton 
et  al  (1971)  and  Dey  et  al  (1975).  This  study,  however,  has 
neglected  any  physical  processes  such  as  precipitation  as 
well  as  orography  and  the  effect  of  an  initial  divergence  in 
these  cases  may  well  be  more  Important. 
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